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Abstract 

The self-consistent Relativistic Quasiparticle Random Phase Approximation (RQRPA) is ex- 
tended by the quasiparticle-phonon couphng (QPC) model using the Quasiparticle Time Blocking 
Approximation (QTBA). The method is formulated in terms of the Bethe-Salpeter equation (BSE) 
in the two-quasiparticle space with an energy-dependent two-quasiparticle residual interaction. 
This equation is solved either in the basis of Dirac states forming the self-consistent solution of 
the ground state or in the momentum representation. Pairing correlations are treated within the 
Bardeen-Cooper-Schrieffer (BCS) model with a monopole-monopole interaction. The same NL3 set 
of the coupling constants generates the Dirac-Hartree-BCS single-quasiparticle spectrum, the static 
part of the residual two-quasiparticle interaction and the quasiparticle-phonon coupling amplitudes. 
A quantitative description of electric dipolc excitations in the chain of tin isotopes (Z = 50) with 
the mass numbers A = 100, 106, 114, 116, 120, and 130 and in the chain of isotones with (A'' = 50) 
88gj.^ ^°Zr, ^^Mo is performed within this framework. 

The RQRPA extended by the coupling to collective vibrations generates spectra with a multitude 
of 2q (g) phonon (two quasiparticles plus phonon) states providing a noticeable fragmentation of the 
giant dipole resonance as well as of the soft dipole mode (pygmy resonance) in the nuclei under 
investigation. The results obtained for the photo absorption cross sections and for the integrated 
contributions of the low-lying strength to the calculated dipole spectra agree very well with the 
available experimental data. 

PACS numbers: 21.10.-k, 21.60.-n, 24.10.Cn, 21.30.Fe, 21.60.Jz, 24.30.Gz 
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I. INTRODUCTION 



Theoretical approaches based on Covariant Density Functional Theory (CDFT) remain un- 
doubtedly among the most successful microscopic descriptions of nuclear structure. The 
CDFT approaches are derived from a Lorentz invariant density functional which connects in 
a consistent way the spin and spatial degrees of freedom in the nucleus. Therefore, it needs 
only a relatively small number of parameters which are adjusted to reproduce a set of bulk 
properties of spherical closed-shell nuclei jj, 121] and it is valid over the entire periodic table. 
Over the years, Relativistic Mean-Field (RMF) models based on the CDFT have been suc- 
cessfully applied to describe ground state properties of finite spherical and deformed nuclei 
over the entire nuclear chart [3] from light nuclei 4] to super-heavy elements j^, ^, from the 
neutron drip line where halo phenomena are observed 7], to the proton drip line with 
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ativistic cranking approxima- 



nuclei unstable against the emission of protons [q]. 
tion has been developed to calculate rotational bands 

excited states, the Relativistic Random Phase Approximation (RRPA) [l^ and the quasi 
particle RRPA (RQRPA) 13| have been formulated as the small amplitude limit of the 



ni l. For a description of nuclear 
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time-dependent RMF models. These models have provided a very good description for the 



positions of giant resonances and a theoretical interpretation of the low-lying dipole [13| and 
quadrupole I3, llSl] excitations. Proton-neutron versions of the RRPA and the RQRPA have 
been developed and successfully applied to the description of spin/isospin excitations as the 
Isobaric Analog Resonance (lAR) or the Gamow- Teller Resonance (GTR) 161 ]. 
Recently, several attempts have been made to extend the RMF and RRPA formalism be- 
yond the mean field approach, first of all, to solve the well known problem of the RMF 
single-particle level density in the vicinity of the Fermi surface which is too low because of 
the too small effective mass. The energy dependence of the single-nucleon self-energy was 
emulated in a phenomenological way 17|] and microscopically by coupling the single particle 
configurations to low- lying surface vibration 1^. This provided a considerable improvement 
for the description of the single-particle spectra. An addition, the quadrupole motion has 
been studied within the relativistic Generator Coordinate Method (GCM) [lol. I20I. 



In Refs. 21 



221 ] ■ we have extended the relativistic RPA by introducing a coupling to col- 



lective vibrations using the techniques developed and realized 



approaches in terms of the Green's function formalism 23 
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ong ago for non-relativistic 



261 ]. An induced additional 



3 



interaction between single-particle and vibrational excitations provided a strong fragmen- 
tation of the pure RRPA states causing the spreading width of giant resonances and the 
redistribution of the pygmy strength to lower energies. This method does not include pair- 
ing correlations and therefore it is restricted essentially to the few nuclei with doubly closed 
shells in the nuclear chart. 

In the present work we consider systems with pairing correlations. Again, we are guided by 
ideas of the quasiparticle time-blocking approximation (QTBA) developed and applied for 
non-relativistic systems in Refs. and 28], which takes into account quasiparticle-phonon 
coupling (QPC) and pairing correlations on an equal footing. However, our approach is based 
on CDFT and formulated in terms of relativistic Green's functions of the Dirac-Hartree- 
Bogoliubov (DHB) or the Dirac-Hartree-BCS (DHBCS) equations. Similar, but in details 
different, app roaches developed earlier within a non-relativistic formalism can be found in 



Refs. 



29 
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3l|. 



32] is that the two types 



The main assumption of the quasiparticle-phonon coupling model 
of elementary excitations - two-quasiparticle and vibrational modes - are coupled in such a 
way that configurations of 2q (8> phonon type with low-lying phonons strongly compete with 
simple 2q configurations close in energy or, in other words, that quasiparticles can emit and 
absorb phonons with rather high probabilities. Obviously, these processes should affect both 
the ground and excited states and therefore, the corresponding amplitudes should be taken 
into account both in the single-nucleon self-energy and in the effective interaction in the 
nuclear interior. 

In order to describe excited states in nuclei, we extend covariant density functional theory by 
coupling the quasiparticles to low-lying vibrations in a consistent way using effective inter- 
actions derived from the same Lagrangian without additional phenomenological parameters. 
First of all, we use the well-known quasiparticle formalism, where, in terms of second quanti- 
zation, nucleon creation and annihilation operators become components of a two-component 
operator mixing a creation and annihilation of a particle into a single quasiparticle. This 
leads to the fact, that for systems with pairing correlations all quantum operators become 
tensors in the two-dimensional quasiparticle space. In particular, the relativistic energy 
unctional is expressed in terms of the relativistic extension of the Valatin density matrix 



33l | of double dimension containing the normal as well as the abnormal densities. As dis- 



cussed in detail in Refs. 



34, 



35l . l36l . l37l | pairing correlations can be considered in a very 
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good approximation as a non-relativistic effect and therefore the full density functional is a 
sum of the relativistic energy functional depending on the normal density and derived from 
the underlying Lagrangian and a non-relativistic pairing energy Epair, depending on the 
abnormal density. The equations of motion are the self-consistent Relativistic Hartree Bo- 
goliubov (RHB) equations. They are derived from this general functional by variation with 
respect to the Valatin density matrix. They are solved numerically and the self-consistent 
fields obtained in this way, which do not depend on the energy, form the static part of the 
nucleon self-energy. This static part determines the nuclear ground state in the mean field 
approximation. 

The static effective interaction used in conventional QRPA approximation is derived as 
the second derivative of the same energy functional and therefore it contains no additional 
parameters. It enables us to go a step further and to compute amplitudes, or vertices, 
which describe the emission or absorption of phonons by quasiparticles within the relativistic 
framework. These amplitudes form the essential ingredient for the following considerations. 
They determine an additive energy- dependent and non-local term in the self-energy of the 
single-quasiparticle equation of motion and, consequently, an induced effective interaction 
between the quasiparticles. Both of these quantities have an influence on the ph- as well as 
on the pp-channel. 

For the calculation of the response of a nucleus in an external fleld we use the Bethe-Salpeter 
equation. It contains both the static and the induced effective interactions and it is for- 
mulated in the doubled two-quasiparticle basis of the Dirac-Hartree-Bogoliubov eigenstates. 
This Bethe-Salpether equation describes the quasiparticle-phonon coupling and pairing corre- 
lations on the equal footing. It is solved using the quasiparticle time blocking approximation 
(QTBA) developed in Ref. 27|], which allows the truncation to 2q ® phonon conflgurations 
and guarantees that the solution is positive deflned. We also use the subtraction procedure 
introduced and justifled in the Ref. j27|. As in the case without pairing it avoids double 
counting of the QPC At zero energy, i.e. at the ground state, particle vibrational coupling 
should have no influence, because the correlations induced by QPC in the ground state have 
already been taken into account in the RHB description through the parameters of the energy 
functional initially fltted to reproduce experimental data, such as nuclear binding energies 
and radii. Therefore, the relativistic mean fleld contains effectively all the correlations in the 
static approximation. The energy dependence of the self energy influences only excitations 
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at finite energy in the nucleus. 

In the present work we develop the Relativistic Quasiparticle Time Blocking Approximation 
(RQTBA) and apply it for the description of electric dipole excitations in even-even spherical 
open-shell nuclei, such as the tin {Z = 50) isotopes ioo,io6,ii4,ii6,i2o,i30gj^ ^j^g ^ 5q) 

isotones ^^Si, ^°Zr, ^^Mo. The RQTBA method, whose physical content is an extension of 
the RQRPA by a coupling to low-lying collective vibrations, provides spectra enriched with 
the 2q ® phonon states. They cause a strong redistribution of the RQRPA strength. As a 
result, we obtain an additional broadening of the giant dipole resonance and a spreading of 
the soft dipole mode (pygmy resonance) to lower energies in the nuclei under investigation. 
The paper is organized as follows. In Section [TTl we formulate basic relations of our approach 
in a rather general form. In Section UTTl we give a more detailed formalism for spherical nuclei 
in the form adopted for numerical calculations. Section [IV] is devoted to the description of 
some numerical details and to the presentation of our results for even-even semi-magic nuclei. 
Finally, Section |V] contains conclusions and an outlook. 

II. GENERAL FORMALISM 

A. Basic relations of the covariant density functional theory for nuclei with pairing 

In this subsection we recall the general formalism of covariant density functional theory with 
pairing, introduce notations and determine conventions used later on. 

In open-shell nuclei, pairing correlations play an essential role and have to be incorporated 
consistently in a description of the ground state as well as of excited states including many- 
body dynamics. Considering pp-correlations in addition to the usual p/i-interaction, existing 
in normal systems, one has to provide a unified description of both pp- and p/i-channels. 
In contrast to Hartree- or Hartree-Fock theory, where pp-correlations are neglected, and 
where the building blocks of excitations (the quasiparticles in the sense of Landau) are 
either nucleons in levels above the Fermi surface (particles) or missing nucleons in levels 
below the Fermi surface (holes), we have now quasiparticles in the sense of Bogoliubov 
which are described by a combination of creation and annihilation operators. This fact can 
be expressed in a standard way by introducing the following two-component operator, which 
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is a generalization of the usual particle annihilation operator: 




(1) 



Here a(l) = a^-^e' is a nucleon annihilation operator in the Heisenberg picture and 
the quantum numbers ki represent an arbitrary basis, 1 = In order to keep the 

notation simple we use in the following 1 = {ri,ti} and omit spin and isospin indices. 
Let us introduce the chronologically ordered product of the operator \E'(1) in Eq. ([1]) and 
its Hermitian conjugated operator ^^(2), averaged over the ground state |$o) of the system 
which will be concretized below. This tensor of rank 2 



G(l,2) = -z(<l>o|r^(l)^1"(2)|$o) 



(2) 



is the generalized Green's function which can be expressed through a 2x2 matrix: 

( a(l)at(2) a(l)a(2) \ 
at(l)at(2) at(l)a(2) 



+ z^(t2-ti)($o| 



at(2)a(l) a{2)a{l] 



(3) 



Similar definitions for t 
used in Refs 



271, 



28 



at(2)at(l) a(2)at(l 
le Green's function in non-relativistic superfluid systems have been 
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391] . Notice that we define the definition of Green's functions here in 



the way of non-relativistic many-body theory, which differs form the conventional definition 
(r\l/^') adopted in relativistic field theories by the replacement of ^ by i.e. by a Dirac 
matrix f3 = 70. This notation is more convenient for our analysis and the matrix (3 needed 
for Lorentz invariance is included in the vertices. Therefore the generalized density matrix 
is obtained ClS db limit 



-i lim ^(1,2) 



(4) 



from the second term of Eq. ([3]), and, in the notation of Valatin 33|, it can be expressed 
as a matrix of doubled dimension containing as components the normal density p and the 
abnormal density x, the so called pairing tensor: 



p{ri,r2,t) x(r-i,r2,t) 
-x*(ri,r2,t) 6{ri - rg) - p*{ri,r2,t) 

These densities play a key role in the description of a superfluid many-body system. 



(5) 



In covariant density functional theory for normal systems the ground state of the nucleus 
is a Slater determinant describing nucleons, which move independently in meson fields (pm 
characterized by their quantum numbers for spin, parity and isospin. In the present investi- 
gation we use the concept of conventional relativistic mean field theory and include the a, u, 
p-meson fields and the electromagnetic field as the minimal set of fields providing a rather 
good quantitative description of bulk and single-particle properties in the nucleus |40|. l4l[|. 
This means that the index m runs over the different types of fields m = {(7,uj, p, A}. The 
summation over m implies in particular scalar products in Minkowski space for the vector 
fields and in isospace for the p-field. In order to obtain a Lorentz invariant theory, these 
classical fields 0^ = {(J,uj'^ , p Af^} are generated in a self-consistent way by the exchange 
of virtual particles, called mesons, and the photon. 

Finally the energy depends in the case without pairing correlations on the normal density 
matrix p and the various fields (prn- 

Ermf[p, <P] = Tr[(ap + l3m)p] + J]|Tr[(/?r^(/)^)p] ± J [^(V^^)^ + U^] d^r^ (6) 

m 

Here we have neglected retardation effects, i.e. time- derivatives of the fields (pm- The plus 
sign in Eq. (jS]) holds for scalar fields and the minus sign for vector fields. The trace operation 
implies a sum over Dirac indices and an integral in coordinate space, a and (3 are Dirac 
matrices and the vertices are given by 

T^ = g^, T^^ = 9^r, f; = g,T^^, F^f = e^^^7^ (7) 

with the corresponding coupling constants gm for the various meson fields and for the elec- 
tromagnetic field. 

The quantities Um{(p) are, in the case of a linear meson couplings, given by the term 



UM) = \ml^l (8) 
containing the meson masses m^- For non-linear meson couplings, as for instance for the 



(T-meson in the parameter set NL3 we have, as proposed in Ref. 42| : 



U[a) = -m^a + —a + —a . (9) 



with two additional coupling constants g2 and g^. 
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In superfluid covariant density functional theory the energy is a functional of the Valatin 
density TZ and the fields (pm- Therefore Relativistic Hartree-Bogoliubov (RHB) theory can be 
derived from an energy functional which depends on the normal density p and the abnormal 
density x as well as on the meson and Coulomb fields (pm- We use here a density functional 
of the form 

ErhbIP, X, X*, 0] = ErmfIP, <P] + Epair[>C, X*] (10) 

where the pairing energy is expressed by an effective interaction V^^ in the pp-channel: 

^pai.[x,x*] = ^rr[x*y^'Px]. (11) 

Here and in the following a tilde sign is used to express the static character of a quantity, 
i.e. the fact that it does not depend on the energy. Of course, in Eq. (fTOj) we could also 
use density dependent pairing forces with Epair = Epair[p,>c] as it is done for instance in 



Refs. 



43 



44| . However, in the present investigation we do not consider this possibility. The 
effective interaction V^^ in the particle-particle channel is supposed to be independent on 
the interaction in the particle-hole channel (see, e.g., Ref. p^]) mediated by the mesons and 
the electromagnetic fields determined above. Generally, the form of V^^ is restricted only 
by the conditions of the relativistic invariance of Epair with respect to the transformations 



of the abnormal densities (see Ref. 



f. 4^). In this section, we consider the general form of 
V^^ as a non-local function in coordinate representation. In all the applications discussed in 
Section IIIII we use for V^^ a simple monopole-monopole interaction. 
The classical variational principle applied to the energy functional of Eq. ( ITOl) 

t2 



6 J (^{<l>o\idt\%) - ErhbIp, X, x*,0])dt = O 



(12) 



leads to the equation of motion for the generalized density matrix TZ: 

tdtn=[nRHB{n),n] (13) 

with the RHB Hamiltonian 

H^.^^.f"^-"-^ ^ ], (14) 

\^ -A* -h^* + m + Xj 

where A is the chemical potential (conted from the continuum limit). In the static case we 
find 

[nRHBin),n] = o. (15) 



Because of time reversal invariance the currents vanish and we obtain the single nucleon 
Dirac Hamiltonian 

= ap + I3{m + S) 

with the RMF self-energy 

S(r) = ^r^0„(r) 

m 

The pairing field A reads in this case: 

A(r,r') = ^ / dr"dr"'VPf{r,r',r",r"')x{r",r"'). 



(16) 
(17) 



Eq. f|T5l) leads to the relativistic Hartree-Bogoliubov equations [35 1 

nnHBm = vEkm, V = ±l 

where 1-0^) are the eigenfunctions corresponding to eigenvalues r]Ek. 
dimensional Bogoliubov-Dirac spinors of the following form 



(19) 

They are the 8- 



l^.^r-)) 



(20) 



Uk{r) 

Vkir) ) \ U*,(r) 

Note that the index k labels here and in the following quasiparticles in contrast to the index 
ki used after Eq. ([T]) for the particle basis. In the following we call this quasiparticle basis 
the Dirac-Hartree-Bogoliubov (DHB) basis. 
The generalized density matrix is obtained as follows: 



7^(r,r') = 5^|^,(r))(^,(0 



(21) 



where the summation is performed only over the states having large upper components of the 
Dirac spinors (i.e. large functions f{k){r) in Eg. f[63l) below). This restriction corresponds 
to the so-called no-sea approximation (see Ref. |37l|). 

The behavior of the meson and Coulomb fields is derived from the energy functional ( fTOi) 
by variation with respect to the fields (pm- We obtain Klein-Gordon equations. In the static 
case they have the form 



AMr) + U'iMr)) = V;(r)/3r^y;(r), 



(22) 



Eq. fl22l) determines the potentials entering the single-nucleon Dirac Hamiltonian f|T6l) and is 
solved self-consistently together with Eq. IHM . The system of Eqs. (HM and (B^ determine 
the ground state of an open-shell nucleus in the relativistic Hartree-Bogoliubov approach. 
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B. Quasiparticle-vibration coupling as a model for an energy dependence of the 
single-quasiparticle self-energy 



The single quasiparticle equation of motion (fT9|) determines the behavior of a nucleon with 
a static self energy. To include dynamics, i.e. a more realistic time dependence in the self 
energy one has to extend the energy functional by an appropriate term leading to a self-energy 
(fTTl) with time dependence. In the present work we use for this purpose the successful but 



relatively simple particle- vibration coupling model introduced in Refs. [32, |46|]. Following the 
general logic of this model, we consider the total single-nucleon self-energy for the Green's 
function defined in Eq. ([2]) as a sum of the RHB self-energy and an energy- dependent 
non-local term in the doubled space: 

E(r-, r'; e) = S(r, r') + J:^'\r, r'; e) (23) 

with 

f (3t(r)6(r -r') A(r,r') \ 
S(r,r') = K ^ ^ ^ ' ^ ' ^ . (24) 

\^ -A*(r,r') -/?S*(r)5(r-r') J 
The energy-dependent operator S*^^^(r, r'; e) will be determined below (the upper index 
e in this quantity indicates the energy dependence). The Dyson equation for the single- 
quasiparticle Green's function ([2]) in the doubled space has the following form: 

{e-n^HB-^^'\^))Gie) = l (25) 

To study the influence of the energy-dependent part of the self-energy on the single- 
quasiparticle energies, it is convenient to formulate Eq. (!25|) in the basis of the eight- 
component Dirac spinors {i/j^) which diagonalize the static RHB-Hamiltonian T-Crhb in Eq. 



?7=±1 k 



where 



Sar^(£) = y rfVdV (V',^Hr)|S(^)(r,r';.)|V^,^^^(r')), (27) 

= J dhd\' «(r)| G{r, r'- eWul{r')). (28) 

In this basis the single-quasiparticle Green's function G of the static mean field has the 
following simple diagonal form: 

GlZ^e) = 5,^,,5,,,M\{el Gl{e) = ^ 5 - +0. (29) 

11 



As in Refs. 

1181, Q, 

we use the particle-phonon coupling model for the energy- dependent 
part of the self-energy Tj^^\ In the basis of the spinors 1-0^) of Eq. ( |20i) . called in the following 
Dirac basis, its matrix elements are given by: 

= E E E '"'f ,n'^'^' ^ ^ - +0- (30) 

The index formally runs over all single-quasiparticle states in the DHB basis including an- 
tiparticle states with negative energies. In the doubled quasiparticle space we can no longer 
distinguish occupied and unoccupied states considering that all the orbits are partially occu- 
pied. But in p ractical calculations, it is assumed that there are no pairing correlations in the 
Dirac sea 33] and the orbits with negative energies are treated in the no-sea approximation. 
As it has been shown in calculations for nuclei with closed shells in Ref. [3], the numerical 
contribution of the diagrams with intermediate states k with negative energy is very small 
due to the large energy denominators in the corresponding terms of the self-energy (|30|) . The 
index /i in Eq. (130|) labels the set of phonons taken into account. are their frequencies 
and 77^ = ±1 labels forward and backward going diagrams in Eq. fl30|) . The vertices '^^J^j^l]^^ 
determine the coupling of the quasiparticles, to the collective state //: 



■ny.mm _ r mm , e .mm* /qi \ 

In the conventional version of the particle- vibrational coupling model the phonon vertices 7^ 
are derived from the corresponding transition densities TZ^ and the static effective interaction: 

kzki rism 

where V^^^^'^^^s denotes a relativistic matrix element of the static residual interaction in 
the doubled space. It is obtained as a functional derivative of the relativistic mean-field 
self-energy S with respect to the relativistic generalized density matrix TZ: 

^'^k2kt 

The transition densities 7^^ are defined by the time dependence of the generalized density 

7^(^) = 7^o + Y^{n^e'''^' + h.c.) (34) 
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describing the oscillating system. We use the linearized version of the model which assumes 
that the transition densities TZ^ are not influenced by the particle-phonon coupling and that 
they can be computed within relativistic QRPA. In the linearized version of the QPC model 
we solve the usual QRPA equations for transition densities 



fca fe4 r)' 

where 

Tp'n _ TpV-V r}(0)v ( \ _ r>{0)v-V ( \ t>W _ -[/V-v' -V,v' (oc\ 

'V;fcifc2 ~ 'V;fcifc2' ^kik2\^) ~ ^kik2 v'^/' ^kikiMk^ ~ ^ kik4,k2k3 ' v""*"/ 

which means that we cut out certain components of the tensors in the quasiparticle space. 
The quantity R is, as usual, the two-quasiparticle propagator, or the mean-field response 
function, which is a convolution of two single-quasiparticle mean-field Green's functions f l29l) : 



'quj — — J^k2 

In Eq. (!35|) we use the static quasiparticle- interaction V of Eq. (!33|) . Of course, in general, 
we should calculate these transition densities taking into account the also the additional 
energy-dependent residual interaction V^^^ [see Eq. (jHj) below] in a self-consistent iteration 
procedure. However, this is not done in the investigations presented here. 



C. Response function in the quasiparticle time-blocking approximation 

Now we have to formulate the Bethe-Salpeter equation (BSE) for the response of a superfluid 
nucleus in a weak external field. The method to derive the BSE for superfluid non-relativistic 



systems from a generating functional is known and can be found, e.g., in Ref. [27| where 
the generalized Green's function formalism was used. Applying the same technique in the 
relativistic case, one obtains a similar ansatz for the BSE. It is formulated now in the basis 
of the DHB spinors in Eq. (!20|) . In full analogy to the case without pairing described in 



Ref. 2l| it is convenient to begin in the time representation. Let us therefore include the 
time variable and the variable t] defined in Eq. (|T9l) . which distinguishes components in the 
doubled quasiparticle space, into the single-quasiparticle indices using 1 = {A;i, t/i, ti}. In 
this notation the BSE for the response function R reads: 

R{U, 23) = G(l, 3)^(4, 2) - * ^ G(l, 5)^(6, 2)y(58, 67)i?(74, 83), (38) 

5678 
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where the summation over the number indices 1,2,... imphes integration over the respective 
time variables. The function G is the exact single-quasiparticle Green's function, and V is 
the amphtude of the effective interaction irreducible in the p/i-channel. This amplitude is 
determined as a variational derivative of the full self-energy S with respect to the exact 
single-quasiparticle Green's function: 



Similar as in Ref. [2l|, we introduce the free response i?°(14,23) = (7(1, 3)G'(4, 2) and 
formulate the Bethe-Salpeter equation ( !38l) in a shorthand notation, omitting the number 
indices: 

R = R^- iR^VR. (40) 

For the sake of simplicity, we will use this shorthand notation in the following discussions. 
Since the self-energy in Eq. fl23l) has two parts S = S + T.^^\ the effective interaction V in 
Eq. (138|) is a sum of the static RMF interaction V and the energy-dependent term V^^^: 

V = V + V^'\ (41) 

where (with ti2 = ti — ^2) 

V{U, 23) = %lZ%7:^(t,,)6{tMt3.) , (42) 

and V^^^'^^l]^ is determined by Eq. (133|) . In the DHB basis of Eq. (l20l) the Fourier transform 
of the amplitude V^^^ has the form: 

W.fc3 ' '^^+0- (44) 



In order to make the Bethe-Salpeter equation (1401) more convenient for the further analysis 
we eliminate the exact Green's function G and rewrite it in terms of the mean field Green's 
function G which is diagonal in the DHB basis. In time representation it has the following 
ansatz: 

(5(1, 2) = -z^^,^,^,^^(r/,r)e-^''^^'=l^ r = h- (45) 
and its Fourier transform is given by Eq. fl29|) . 
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Using the connection between the mean field GF G and the exact GF G in the Nambu form 

G-i(l,2) = G-i(l,2) + S^(l,2), (46) 

one can ehminate the unknown exact GF G from the Eq. (HOj) and rewrite it as follows: 

R = R^- iR^WR (47) 

with the mean-field response i?°(14, 23) = G{1, 3)G'(4, 2), and is a new interaction of the 
form 

W = V + W''^\ (48) 

where 

W^'\U, 23) = V^'\U, 23) + 2S(^)(1, 3)G-^(4, 2) + iG~\l, 3)S(")(4, 2) - 3)J:^'\A, 2). 

(49) 

Thus, we have obtained the BSE in terms of the mean-field propagator, containing the well- 
known mean-field Green's functions G, and a rather complicated effective interaction W in 
Eq. (H8|) . which, however, is also expressed through the mean- field Green's functions. 
The structure of the energy-dependent effective interaction W^^^ has a clear interpretation 
in terms of Feynman's diagrams which are usually employed to clarify the physical content 
of the amplitude Vl^^^^ {2^, . In addition to the static interaction V, the effective inter- 
action W contains diagrams with energy- dependent self-energies and an energy- dependent 
induced interaction, where a phonon is exchanged between the two quasiparticles. In the 



present work, as well as in Ref. 21J, we omit the term zS(^^(3, 1)S'^'^)(2, 4) in Eq. fH9l) 
because it plays a compensational role with respect to the backward-going components of 
the previous terms in the W^'^\ However, within the version of the time blocking ap- 
proximation, which we apply to the BSE (see below), the backward-going propagators are 
not taken into account. Components containing the backward-going propagators within 
2q ® phonon configurations require a special consideration which is formulated in Ref. 
for a superfiuid non-relativistic system. In the present work these correlations are fully ne- 
glected and, therefore, the term zS^'^^(3, 1)S'^^)(2, 4) has also to be omitted. However, we 
have to emphasize, that we only neglect ground state correlations (GSC) (backward-going 
diagrams) caused by the quasiparticle-phonon coupling. All the QRPA ground state cor- 
relations are taken into account, because it is well known that they play a central role for 
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the conservation of currents and sum rules. We consider that this is a reasonable approx- 



imation w 
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lich is a pp. 



27 



28 



31 



led and discussed also in some non-relativistic models (see e.g. Refs. 



47 



49 



50 



5l| and references therein). 



Eq. P7j) whose integral part contains singularities in the amplitude W can not be solved 
explicitly because, considering the Fourier transform of the Eq. ( 147|) . one finds that both the 
solution of this equation R and its kernel W are singular with respect to energy variables. 
Also, this equation contains integrations over all time points of the intermediate states. This 
implies that many configurations which are actually more complex than 2q ® phonon are 
contained in the exact response function. Therefore, we apply the special time-projection 
technique, introduced in the Ref. 2J] and generahzed in Ref. |27(] for superfluid systems, to 
block the 2q'-propagation through these complicated intermediate states. 
Conventionally, we divide the problem to find the exact response function of the BSE (1471) 
into two parts. First, we calculate the correlated propagator R^^^ which describes the 2q- 
propagation under the influence of the interaction 



r(^) = R^ - iR^W^^^ R^^\ 



(50) 



It contains all the effects of particle-phonon coupling and all the singularities of the integral 
part of the initial BSE. Second, we have to solve the remaining equation for the full response 
function R 

R = R(^^ -iR^^^VR. (51) 

Eq. (EH) contains only the static effective interaction V and can be easily solved when R^^^ 
is known. 

The correlated propagator R^'^^ can be represented as an infinite series of graphs which contain 
mean-field 2g-propagators alternated with single interaction acts. This can be expressed by 
the system of the following equations employing the auxiliary amplitude F^*^^: 



i^(^) = /^o-^/^°F(^)/^^ 



(52) 
(53) 



Then, the integral part of Eq. (1531) has to be modified to order in time the interaction acts 
described by the amplitude W^^\ It means that we should cut out only terms where the 
'left' time arguments of the amplitude F^*^^ are greater than the 'right' time arguments of the 
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amplitude IV^^^ |2J,l27|. This can be expressed by the time-projection operator of the form: 

6(14,23) = Sn^_r^2^kik3Sk2kiO{r]itu)0{r]it23), (54) 
which is introduced into the integral part of the Eq. fl53|) : 

r('^)(14, 23) = M^(")(14, 23) + - ^ W^'^^IG, 25)^°(58, 67)9(58, 67)r(")(74, 83). (55) 



5678 



Since we are interested in spectral characteristics of the nuclear response, a Fourier trans- 
formation of the response function is performed as follows: 

oo 

KZmI>) = j dhdt2dhdU6{h - t2)5{h - t4)5(t4)e*"*"/2(14, 23), (56) 

— oo 

so that the response function depends only on one energy variable uj. 

The time projection by the operator ( 154|) leads, after some algebra and the transformation 
(1561) . to an algebraic equation for the response function. For the p/i-type components of the 
response function it has the form: 

where 

^kik4,k2k3i^) ~ ^kik4,k2k3 ik4,k2k3 

{0))S,,,. (58) 

$ is the particle-phonon coupling amplitude in the QTBA with the following forward (r/ = 1) 
and backward {i] = —1) components: 



-Tj —Tj* Tj Tj* 



ka ' p 

—n 

(59) 



IJ,;kik3^ ^]k2k4 ^ T/x;fc3A;iT/x;fc4fc2 

— — Ek^ — r]u — Ek^ — Ek^ — il^, . 

where we denote: 

^1-Mk2 = ^]^-klk2 ^tk4,k2ks('^) = ^t~klk2Vi^)- (60) 

Indices ki in this expression formally run over the whole DHB space, but in applications 
we usually consider that the amplitude ^2ik4 k2kS'^) describes phonon coupling only within 
some energy window around the Fermi surface. That is why it implies that this amplitude 
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contains no antiparticle-quasiparticle (aq) configurations. Notice, that in our approach we 
cut out only the components without ground state correlations induced by phonon coupling 
(l60ll which include the main contribution of the phonon coupling and neglect some more 
delicate terms. 

However, ground state correlations of the QRPA type are taken into account due to the 
presence of the V^^f^ ^^ks terms of the static interaction in the Eq. flFTI) . By definition, the 
propagator R{uj) in Eq. (1571) contains only configurations which are not more complicated 
than 2q (g> phonon. 

In Eq. (1571) w e have included the subtraction procedure because of the same reasons as 
in the Ref. 2l|. Since the RMF ground state is adjusted to experimental data, it contains 
effectively many correlations in the static approximation and, in particular, also admixtures 
of phonons. Therefore, when we include them explicitly in the dynamics, this static part 
should be subtracted from the effective interaction to avoid double counting of the QPC 
correlations. Since the parameters of the density functional and, as a consequence, the 
effective interaction V are adjusted to experimental ground state properties at the energy 
00 = 0, this part of the interaction ^{uj), which is already contained in V, is given by $(0). 
This subtraction method has been introduced in the Ref. [27| for self-consistent schemes. 
Eventually, to describe the observed spectrum of the excited nucleus in a weak external field 
P as, for instance, an electromagnetic field, one needs to calculate the strength function: 

S(E) = -- lim ImU(E + iA), (61) 

expressed through the polarizability n(co') defined as 

The imaginary part A of the energy variable is introduced for convenience in order to obtain 
a more smoothed envelope of the spectrum. This parameter has the meaning of an addi- 
tional artificial width for each excitation. This width emulates effectively contributions from 
configurations which are not taken into account explicitly in our approach. 
In relativistic RPA and QRPA calculations the Dirac sea plays an important role. A consis- 
tent derivation of relativistic RPA (QRPA) as the small amplitude limit of time-dependent 
RMF (RHB) theory in Ref. |l^ shows that one has to include besides the usual ph- 
configurations also antiparticle-hole (ah) configurations. Otherwise current conservation 
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is violated 5^ and the position of giant resonances cannot be described properly in rel- 



ativistic RPA 53|]. However, this increases the number of configurations dramatically as 



compared to non-relativistic QRPA calculations and requires in particular in deformed rel- 



ativistic QRPA calculations 



54j a tremendous large numerical effort. Recently a simple 



method has been proposed to avoid this problem. As discussed in Ref. [55|], the static no- 
sea (SNS) approximation takes the contributions of the empty Dirac sea into account in a 
very good approximation by a renormalization of the total effective interaction W{lj) in the 
Bethe-Salpeter equation. 

III. APPLICATION OF THE APPROACH: BASIC APPROXIMATIONS 

The formulated relativistic QTBA is applied to calculations of the dipole strength in spherical 
nuclei with pairing. In this application we mainly follow the calculation scheme employed 



in Ref. 
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22| . however, with some considerable modifications accounting pairing effects: 
all the equations are solved in the doubled space. The computation is performed by the 
following main steps: 

i) To calculate ground state properties the Dirac equation together with the BCS equation 
for single nucleons are solved simultaneously with the Klein-Gordon equation for meson fields 
in a self-consistent way to obtain the single-quasiparticle basis (Dirac-Hartree-BCS basis). 

ii) The RQRPA equations (l35l) with the static interaction V of Eq. (1331) are solved in the 
Dirac-Hartree-BCS basis to determine the low- lying collective vibrations (phonons), their 
energies and amplitudes. In the present work we have included the phonon modes with 
energies below the neutron separation energies for the Z=50 chain and with energies below 
10 MeV for the N=50 chain. The two sets of quasiparticles and phonons form the multitude 
of 2q® phonon configurations which enter the quasiparticle-phonon coupling amplitude $(a;) 
in Eq. (|59l). 

iii) The equation for the response function ( l47l) is solved using this additional amplitude in the 
effective interaction W{(jj). Making a double convolution of the response function with the 
external field operator P, one obtains the polarizability (l62ll and the strength function ( 16T|1 
determining the spectrum of the nucleus. It is found that the amplitude W{uj), containing 
a large number of poles of 2q ® phonon nature, provides a considerable enrichment of the 
calculated spectrum as compared to the pure RQRPA. 
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A. Description of the ground state 



In the present work we confine ourselves by the case of spherically symmetric nuclei where it is 
convenient to separate the dependence on the magnetic quantum number ruk'- k = {{k), ruk}, 
where (k) is the set of remaining quantum numbers which are time reversal invariant: —k = 
{{k), —nik}- In this case {k) = {rik, jk,nk,Tk} with the radial quantum number n^, angular 
momentum quantum number j^, parity tt^ and isospin r^, so the Dirac spinors read: 

f ,N I f{k)ir)yi,j,m,i^,ip) \ 
Mr,t) = \ \Xr,{t), (63) 

yijm{'>^, v) is a two-component spinor 

3^0m(^' V^' = ^C-^slm\j'rn)Yi^X'^,ip)xrns{s) , (64) 

t is the coordinate for the isospin and Xr^i'l^) is a spinor in the isospin space. The orbital 
angular momenta Ik and Ik of the large and small components are determined by the parity 
of the state k: 

k = 3k + |, 4 = jk-\ for TTfc = (-l)J'=+5 

1 (oSj 
Ik = jk - |, 4 = jk + \ for VTfc = (-l)-"'-2, 

/(fc)(r) and g{k){f) are radial wave functions. The phase convention for the wave function 
is chosen so that the following relation is fulfilled: 

l^l\l = {-)'''^'''-'^^^-k- (66) 

In the literature [l^ the RQRPA are solved for finite range Gogny forces in the pairing chan- 
nel in the canonical basis. This has the advantage, that the quasiparticle matrix elements of 
the QRPA-equations can be calculated rather easily by multiplying the matrix elements in 
particle space by BCS-occupation factors, but it has the disadvantage, that the matrix H^^ 
in quasiparticle space is no longer diagonal in the canonical basis. The quasiparticle energies 

+ have to be replaced by complicated matrices. 
We therefore use in the following applications the RMF-f-BCS approximation, where the 
canonical basis coincides with the BCS-basis. In this approximation the ground state wave 
function |$o) is considered to be a vacuum state with respect to quasiparticles with the 
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creation and annihilation operators a^, at determined by the special Bogoliubov transfor- 



mat ion: 

ak \ _ I Uk -Vk 
a\ j \Vk Uk 




ak\%) = Vfc, (67) 



where u\ + vf, = 1. Operation k transforms the state k to the time reversal state. In a 
spherical system we define 

a-k = i-iy^+^^-'^''a_k, (68) 
where the choice of the phase factors is determined by Eq. (166|) . 

In the RMF+BCS approximation we determine, in each step of the iteration, first the eigen 
functions (fk of the single-particle Dirac Hamiltonian of Eq. (fT6i) 



dx'h {x,x')(fk{x') = {m + ek)'^k{x) (69) 

where the coordinate x = {r,a,t} combines the spatial coordinates r with the Dirac index 
a = 1 ... 4 and the isospin t. Next the Dirac spinors ipk are used to construct the single- 
particle density matrix 

p{x, x') = '^<fk{x)vl<fl{x'). (70) 

k 

In the basis of the functions (fk (BCS basis) p as well as are diagonal with the eigenvalues 
vl and m + 6k- The pairing field A is in this basis close to canonical form: Akki = Skk'^k- All 
other matrix elements vanish in the case of a monopole force with constant matrix elements 
and without cut-off, in other cases they are neglected in the BCS approximation. Thus, in 
this basis, the Hartree-Bogoliubov matrix (fT9l) is reduced to a set of 2x2 matrices, which can 
be diagonalized analytically. Thus one finds as eigenvalues the quasiparticle energies 



Ek = \Kek-K? + /^l (71) 



and as eigen functions the occupation amplitudes Uk and Vk with 

1 A Sk- Kk 



and Uk = \/ 1 — v\. The pairing gaps are obtained by the solution of the gap equation 

Afc/ 



^k = --yv^i^,r,^ (73) 

k' " 
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in each step of the iteration and the chemical potential A^^^. is fixed via particle number 
conservation: 

vl = N (or Z) for neutrons (or protons). (74) 

k 

After the solution of the BCS equations fl73ti74p the density (1701) is calculated and used for 
the solution of the Klein-Gordon equations fl22|) determining the RMF potentials for the 
Dirac-Hartree Hamiltonian in Eq. fl69|) in the next step of the iteration. In the RMF+BCS 
approximation the eight components of the quasiparticle eight-component Dirac spinor lip^) 
are simply expressed through the usual 4-component spinor wave functions ipk- 

Vk{x) = (-l)''=+^''=+'"'=^;fc(^l,(x), (75) 

and we have chosen Uk,Vk > V/c. This simplifies the calculation of the quasi-particle RPA 
matrix elements in the next section considerably. We only have to calculate the matrix 
elements in particle space using the wave functions (fk and multiply them with the corre- 
sponding BCS occupation factors in Eqs. flHT]) and fl82|) . 

In the present applications of our approach we use a monopole force with constant matrix 
elements and a soft pairing window. Details are given below. 



B. Solution of the RQRPA equations and calculation of the phonon vertices 



The RQRPA equations are derived as the small amplitude limit of the time-dependent Dirac- 
Hartree-Bogoliubov equations for the generalized density matrix TZ ISj. For general pairing 
forces, as for instance for the finite range Gogny force in the pairing channel 561] they 
can be solved in the canonical basis 57| of the RHB equations, where the full Hartree- 
Bogoliubov ground state wave function has BCS form. In the RMR-I-BCS case they are 
solved in the Dirac-Hartree-BCS basis (!75|) described above. In spherical systems we can 
use angular momentum coupling of the 2-quasiparticle states and the reduced form of the 
RQRPA equation for angular momentum is: 

(^^M ~ - Ek2)T^l(^kik2) ^ ^2 ^2 ^(kiki,k2k3)'^l(k3k4)' C^^) 

V' (fc4)<(fc3) 

where the index // characterizes the various solutions of the RQRPA equation, in particular 
their angular momentum J^. The notation {kik2) indicates the fact that the two quasiparti- 
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cles with the indices ki and k,2 are coupled to angular momentum J^. The i]!]' components 
of the static residual interaction in the p/i-channel read: 



S=0,1 



X 



S S ~{ph)JS I tS ~{pp)J 

lfc2,fc3fc4) 



(77) 



where '^(fcl'fc^ fj^a) ^^"^ ^(^^1:2 ^3^4) reduced matrix elements of the ph- and pp-interaction. 

We assume that the pp-components do not depend on the total spin, and the p/i-components 
carry spin S = 0, 1. The p/i-components "^(1423)^ describe the one-boson exchange (OBE) 

interaction and could be expressed as follows: 

00 



(fcifc4,fc2fc3) 2J -\- 

~ m£S L 



lllllj ^^^((^i)ll^-a^0[/5r..F.]^||(A:,))x 

X Di{q,q'){{k,)\\3L{q'rmr.sYLYm)), (7J 



where in the first sum (m G S) the index m runs over the various meson fields carrying 
spin 5". The index S in Tms denotes the spin of the Pauli matrix entering the vertices Tm 
in Eq. ([71). This implies in particular that 5 = for the scalar and time- like parts of the 
vector mesons and that = 1 for the space-like parts of the vector mesons (current-current 
interactions) . 

Representing the g-integral in Eq. flTHl) by a discrete sum over mesh points, the matrix 
elements (|78l) are a sum of separable terms. The non-local meson propagator is a solution 
of the integral equation: 

q'Di{q,q') + / 1^ Mi{q - q")Di{q",q') = {2nf6{q - q'), (79) 

where M^{q) is the Fourier transform of U" {(j)^{r)) determined by Eq. (l8|9l) : 

M^{q) = jd\e-^'^^U"{<Pi{r)). (80) 

The quantities ^(^^^2)' ^(^11^2) -^'l- conventional factors |57i] which are the 

following linear combinations of the occupation numbers: 

' (81) 



^ffcifc2) = /I \ (uk.Uk^ - i~lfvk,Vk2), 



(82) 
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arising due to symmetrization in the integral part of the Eq. (jTHD, which enables one to 
take each 2g-pair into account only once because of the symmetry properties of the reduced 
matrix elements '^^(fc^'fc^fcgfca) ''^(fcffe ^3^4) ' "^^^ interaction v^^^'^ in theop-channel we use 
a simple monopole-monopole ansatz with the so-called smooth window 58| : 



where is the value of the pairing window and d is its diffuseness. 

The RQRPA transition densities ^^^(^j^;,^) calculated from the Eq. (!76l) determine in particu- 
lar the components of the amplitudes 7^(fc^fc2) ^^i^h couple the phonon with the quasiparticle 
states {ipkl) ^^cl having rji = 7]2 = r] i.e. lying on the same side with respect to the 

Fermi level 



ri' (fe4)<{fc3) 5=0,1 



X 



\tS „5 ~(ph)Jt.S _ S pS ~{PP)J,. l-i-,r?' (n.^ 
[?(A;ifc2)'/(fc3fc4)''(fcifc4,fe2A:3) ''(fciA:2)?(fc3fc4)'^(fcife2,A:3fc4)J '^fiik-jki)- 



C. The RQTBA correlated propagator and the strength function 

In solving Eq. (H71) for the response function, we use our previous experience with calculations 



for nuclei with closed shells 
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22| . Again, we formulate and solve this equation both in the 
2g-basis of Dirac-Hartree-BCS quasiparticle pairs and in the momemtum-channel space. In 
Dirac-Hartree-BCS space its dimension is the number of 2g-pairs which satisfy the selection 
rules for the given multipolarity. In relativistic nuclear calculations it is always important 
to take into account the contribution of the Dirac sea. This can be done, as it is done 



traditionally, exp 
proposed in Ref. 



icitly, or statically by the renormalization of the static interaction, as it is 



551 ]. Nevertheless, for systems with pairing correlations the total number of 



2g-pairs entering Eq. ( H71) increases considerably not only with the nuclear mass number, but 



also with the pairing window. As it was investigated in a series of RRPA calculations 12, 1591], 
the completeness of the ph (ah) basis is very important for calculations of giant resonance 
characteristics as well as for current conservation and a proper treatment of sjTumetries, in 
particular, the dipole spurious state originating from the violation of translation symmetry 
on the mean field level. On the other hand, the use of a large basis requires a considerable 
numerical effort and, therefore, it is reasonable to solve the Eq. ( 1471) in a different more 
appropriate representation. 
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Our choice is determined by the following properties of the static effective interaction V. 
Its p/i-component is based on the exchange of mesons and explicitly contains only direct 
terms and no exchange terms, therefore it can be written as a sum of separable interactions 
(1781) . and in the present work its jop-component is also chosen in the separable form fl83l) for 
convenience. 
As in Refs. 
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22l |. we solve the response equation for a fixed value of the energy variable 



uj in two steps. First, we calculate the correlated propagator R'^{uj) which describes the 
propagation under the influence of the interaction $(a') in the time-blocking approximation 
without GSC caused by the phonon coupling: 



/ ^ [ {kike, 



(uj) — <|)^*^"^''' 



(0) 



r>ie)J,ri / \ 

^ik5k4,kek3)y^)^ 



^5) 



(fc6<fe5) 



where the symmetrized matrix elements of the mean field propagator R^^^ and the two 
quasiparticles-phonon coupling amplitude ^^'^^ read: 



^tklklk2k3)(^) = ^Sl2)(^)('^(felfe3)^(fc2fc4) + (-)'^^'' ^'<^(fclfc4)^(fe2fc3))> 



^ {kik4,,k2k:j) ('^) 



1 



^ (kikiMks)^^ ^ ^ ^ > {k2kiMkay^ 



(86) 
(87) 



1 + ^(kski) 

which means that we take into account two kinds of components: one kind with only forward 
{rj = 1) 2g-propagators of the ph-tjpe {rji = —7^2) and another one with only backward 
propagators (77 = —1), but do not include mixed ones. In the conventional terminology it 
means that we neglect ground state correlations caused by the quasiparticle-phonon coupling. 



The reduced matrix elements of the quasiparticle-phonon coupling amplitude $('^'1^^ 



k2k3) 



read: 



^ {kik4,k2k3)^^^ 



E 



C,k2)'^ ti{kfik4) 



"^^^2 + 1 TZiV^- Eki- Eke 



Tj rj* 
^ ti{kiks)^ fj.{k3k5) 

IV. - Ek2 - fi/. 



+ (-1 



Jk, jk2 J 1 / i-iy'' '''^l{k,k3)^J!k2k4) 
Jk4 Jk3 J,}^ - - Ek2 - % 



+ 



r]uj - Ek^ - E, 



ki 



(88) 
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where denotes the relativistic quantum number set: Xjt = {2jk + ^){h~ jk)- The reduced 
matrix elements of the particle-phonon coupling amplitude T^i^^j^^) calculated from the 
Eq. (IH^ . The index /i = {J^,n^} denotes the set of phonon quantum numbers which are its 
angular momentum and the number of the solution of the Eq. fl76l) . The quantity 
is the corresponding energy. The fact, that the r.h.s. of Eq. (188!) depends only on the same 
r^-values as the l.h.s. and does not contain any mixing of different //-values implies that no 
GSC are contained in the intermediate 2q ® phonon propagators. 

The Eq. ( I85l) is too expensive numerically to be solved in the full Dirac-Hartree-BCS ba- 
sis. However, due to the pole structure of the ^-amplitude it is naturally to suggest that 
quasiparticle-phonon coupling effects are not important quantitatively far from the Fermi 
surface. In the present work, for numerical calculations an energy window Emn was imple- 
mented around the Fermi surface with respect to pure two-quasiparticle energies i?2g so that 
the summation in the Eq. ( l85l) is performed only among the 2g-pairs with £'213 < Emn- Con- 
sequently, the correlated propagator differs from the mean field propagator only within this 
window. This approximation has been checked in the Ref. [21l] in the calculations for nuclei 
with closed shells by direct calculations with different values of this energy window, and it 
has been found that this window should include just the investigated energy region. Beyond 
the energy window we do not obtain additional poles caused by 2q ® phonon configurations, 
but only the renormalized QRPA spectrum. It is important to emphasize that many 2q- and 
aq'-configurations outside of the window are taken into account on the RQRPA level that 
is necessary in order to obtain the reasonable centroid positions of giant resonances as well 
as to find the dipole spurious state close to zero energy. By its physical meaning, the Eq. 
( !85l) contains all effects of the quasiparticle-phonon coupling and all the singularities of the 
integral part of the initial BSE. 

In the second step, we have to solve the remaining equation for the full response function 

R{uy. 

(fc6<fc5) {ke,<kj)r)" 

In contrast to the Eq. fl85|) . this equation contains only the static effective interaction V from 
the Eq. ([771). 

Since both the one-boson exchange interaction and the pairing interaction are separable 
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in momentum space, we can use this advantage and formulate the response equation in 
the momentum-channel representation. Let us introduce the following generalized channel 
index = {q,m, L, S} for x = (ph) and = for x = (pp). For x = ijph) it includes 
the momentum q transferred in the exchange process of the corresponding meson labeled 
by the index m. The index x distinguishes ph- and pp-channel components of the static 
interaction, L is the angular momentum, and the index S=0,1 has its usual meaning of the 
total spin carried through the certain channel. In this way, we apply the following ansatz 
for the T^-components of the static effective interaction V 

cc' 

where we omit the index x for simplicity. For the channels with x = (ph): 

QtZ = ^'''l^^^~' vfk.k.)mn^^^^ (91) 

dec' — ± „ T , 1 — \fi — OLL'OsS'Omm' [^^) 
ZJ -\- 1 (ZTTj" 

and the summation over c, c' implies integration over d^q, d^q' . For the channels with x = (pp) 
we have: 



Qik^k.) - ^Jodik.>^.) ^/TTIST ^'^'^^ V 1 + e(^^i-)/'^ ^^^^ 
G 

dec' = --Scc'- (94) 

Then, we can use the well known techniques of the response formalism with separable inter- 
actions (see, for instance, Ref. [s^])- We define the exact response function and the correlated 
propagator in the generalized momentum-channel space as follows: 

Ri'H= E E QS:S^S1....3)(-)Q(:i') (95) 

(k2<kiyn {ki<ks)-n' 

{A:2<fci) (fc4<A:3)»? 

In this representation Eq. (l89ll reads: 



R^^, = Rl^, + {R'dR)ec'. (97) 
This equation is solved by matrix inversion 

R= (1- R'^d^ ^R\ (98) 
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To compute the nuclear response in the certain external field, we need a convolution of 
the exact response function with the external field operator P which can be suggested as 
an additional channel c = p, p = {z,x}, where the index z contains possible additional 
dependences of the external field which we do not concretize here: 

„x, _ J2 ''■"'Vrr'""' 4.fa)((*0 II II (fe)>- (99) 
Making use of this definition, we can determine the polarizability as: 

n^(a;) = R'juj) = Rlfiu) + ^ i?g)^(^)4c'i?i(a;), (100) 

cc' 

where the quantities E^c\uj), E^p\lo) can be found as follows: 

TD{e)j ( \ _ \^ \^ p{-p)J:n* -pi'^)J,v (, An'^^')-^'^ 

(fc2<fci) (fc4<fc3)r? 

4f(-)= E E (101) 

(fc2<fcl) {k4<ks)ri 

and the quantity R'^p{u), which has a meaning of the density matrix variation in the external 
field P, obeys the equation: 

RiM = RS'(^) + E RlS'i^)dc'c"RUu^). (102) 

c'c" 

To describe the observed spectrum of the excited nucleus in a weak external field P, as for 
instance a dipole field, one needs to calculate the strength function: 

S'^(E) = -- lim ImU'^iE + iA), (103) 

expressed through the polarizability Il'^{u) defined by Eq. fllOOl) . 

Obviously, the dimension of vectors and matrices entering Eq. (11021) is determined by the 
number of mesh-points in g-space and the number of m, L, iS-channels. In particular, it 
does not depend considerably on the total dimension of 2q- and ag-subspaces and on the 
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221, the 



mass number of the nucleus. As we have realized in the calculations of Refs. 
advantage of the momentum- channel representation appears at some medium values of the 
nuclear mass number, where the total dimension of ph- and a/i-subspaces, which is exactly 
the dimension of arrays in the Eq. (1571) written in the coupled form, become comparable with 
the dimension of matrices entering Eq. f ll02p . In the present approach, due to the pairing 
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correlations, this mass region shifts towards lower masses. The solution in the momentum- 
channel space is even more helpful when we include pairing correlations, since the number 
of states within the pairing window increases with more than a factor two as compared to 
the case without pairing. For heavy nuclei the dimension of the two-quasiparticle DHB basis 
increases considerably and, therefore, for heavy nuclei the solution of the response equations 
in momentum space is recommendable. 

Notice, that the pairing correlations cause also an additional numerical effort in Eq. (I85p . 
It is solved within the subspace of 2q'-configurations confined by the which, in the 

realistic calculations, surrounds the pairing window and, therefore, contains considerably 
more configurations as compared to the case of no pairing. 



IV. COMPUTATIONAL DETAILS, RESULTS AND DISCUSSION 



A. Numerical details 



For this first application we have chosen two chains of spherical even-even semi-magic nuclei: 
one chain with Z = 50 and another one with = 50. We have calculated the isovector 
dipole spectrum in the giant dipole resonance region and in the low-lying energy region in 
the two approximations: RQRPA and RQTBA for the quasiparticle-vibration coupling. All 
.he results presented below have been obtained with making use of the NL3 parameter set 



60l | for the covariant density functional (jH]). 

In the present work, pairing correlations were treated in the BCS approximation where 
the single quasiparticle wave functions diagonalize the single-nucleon density matrix p. As 
pairing interaction V^^ we use the simple monopole-monopole form (183!) within the smoothed 
energy window with the parameters w = 20 MeV, d = 1 MeV. The parameter G was chosen 
in such a way that the resulting gap at the Fermi surface reproduces the empirical gap 
expressed by the well known three-point formula: 

= -^-^[B{N, - 1) + B{K + 1) - 2B{Nr)], (104) 

where B{Nt-) is the experimentally known binding energy of the nucleus with A^,- nucleons in 
the subsystem with pairing correlations (neutrons or protons). The RMF plus BCS equations 
are solved by expanding the nucleon spinors in a spherical harmonic oscillator basis [3|]. In 
the present calculation we have used the basis of 20 oscillator shells. 
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In solving the RQRPA Eq. (1751) we have used the method proposed in Ref. 6l| for a reduction 
of the eigenvalue problem by the generalized Cholesky decomposition. In the RQRPA as 
well as RQTBA calculations both Fermi and Dirac subspaces were truncated at energies far 
away from the Fermi surface: in the present work as well as in the Refs. 2l|, we fix the 
limits < 100 MeV and Eaq > —1800 MeV with respect to the positive continuum (so far 
from the Fermi surface there are no pairing effects, therefore we have there pure particles and 
holes). A small artificial width was introduced as an imaginary part of the energy variable 
huj to have a smooth envelope of the calculated curves. In the calculations for tin isotopes 
we took 200 keV smearing for the spectrum in the wide energy region 0-30 MeV and 20 
keV for the low-lying portion of the same spectrum below 10 MeV to distinguish its fine 
structure. For the N=50 isotopes we have used the smearing 400 keV, assuming the more 
pronounced contribution of the single-particle continuum in the GDR region, and 10 keV for 
the low-lying strength. 

The energies and amplitudes of the most collective phonon modes with spin and parity 



2"*", 3 , 4+, 5 ,6"*" have been ca^ 
same criterion as in the Ref. 21 



culated with the same restrictions and selected using the 



22l | and in many other non-relativistic investigations in 



this context. Only the phonons with energies below the neutron separation energy for the 
examined tin isotopes and below 10 MeV - for the = 50 nuclei enter the phonon space 
since the contributions of the higher-lying modes are su ppo sed to be small. Our previous 



experience within the non-relativistic approach of Ref. [28| without the restriction of the 
phonon space by the energy have shown that the inclusion of the high-lying modes into 
the phonon space cause the change of the mean energies and widths of the resonances 
comparable with the smearing parameter (imaginary part of the energy variable) used in 
the calculations, because the physical sense of this parameter is to emulate contributions of 
remaining configurations which are not taken into account explicitly. 

As a test of numerical correctness of our codes, the response equation has been solved both in 
the DHBCS basis and in momentum- channel space and identical results have been obtained. 
Since the quasiparticle-phonon couphng amplitude fl88|) has a pole structure, its contributions 
to the final result for the strength function decrease considerably when we go away from the 
Fermi surface. Therefore, this coupling has been taken into account only within the 2q- 
energy window < 25 MeV around the Fermi surface. This restriction means that above 
this energy we have no poles induced by the complex configurations, and obtain the pure 
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RQRPA poles, but with larger strength which comes from the integral contribution of the 
lower-lying energy spectrum. It has been checked that a further increase of this window does 
not influence considerably the strength functions at energies below the value of this window. 
Although a large number of configurations of the 2q C?> phonon type are taken into account 
exphcitly in our approach, nevertheless we stay in the same two-quasiparticle space as in 
the RQRPA, therefore the problem of completeness of the phonon basis does not arise and, 
therefore, the phonon subspace and the subspace of the 2q® phonon states can be truncated 
in the above mentioned way. Another essential point is, that on all three stages of our 
calculations the same relativistic nucleon-nucleon static interaction V has been employed. 
The vertices (18^ entering the QPC energy-dependent interaction are calculated with the 
same force. Therefore no further parameters are needed, and our calculation scheme is fully 
consistent. 

The subtraction procedure developed in the Ref. j27| for self-consistent schemes has been 
incorporated in our approach. As it was mentioned above, this procedure removes the 
static contribution of the quasiparticle-phonon coupling from the static interaction in the 
p/i-channel. Therefore, the QPC interaction takes into account only the additional energy 
dependence introduced by the dynamics of the system. It has been found in the present 
calculations as well as in the calculations of the Ref. |28| that within the relatively large 
energy interval (0 - 30 MeV) the subtraction procedure provides a rather small increase of 
the mean energy of the giant dipole resonance (about 0.7 MeV for tin region) and gives 
rise to the change by a few percents in the sum rule. This procedure restores the response 
at zero energy and, therefore, it does not disturb the symmetry properties of the RQRPA 
calculations. The zero energy modes connected with the spontaneous symmetry breaking in 
the mean field solutions, as, for instance, the translational mode in the dipole case, remain 
at exactly the same positions after the inclusion of the quasiparticle-vibration coupling. In 
practice, however, because of the limited number of oscillator shells in our calculations this 
state is found already in the RQRPA without the QPC at a few hundreds keV above zero. 
In cases, where the results depend strongly on a proper separation of this spurious state, 
as, for instance, for investigations of the pygmy dipole resonance in neutron rich systems, 
we have to include a large number of the 2g-configurations in the RQRPA solution to avoid 
mixing of the spurious state with the low-lying physical states. 
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B. Isovector dipole strength distribution in semi-magic nuclei: pygmy and giant 
resonances 



In Figs. [T] and [2] the calculated dipole spectra for the tin isotopes ^°°Sn, ^°^Sn, ^^^Sn and 
^^^Sn, ^^°Sn, ^^°Sn, respectively, are given. The right panels of the figures show the photo 
absorption cross section 

^Ei{E) = Sei{E), (105) 

which is determined by the dipole strength function Sei, calculated with the usual isovector 
dipole operator. The left panels show the low-lying parts of the corresponding spectrum 
in terms of the strength function, calculated with the small imaginary part for the energy 
variable, in order to see the fine structure of the spectrum and sometimes individual levels 
in this region. Fig. [3] represents the analogous results for the three = 50 nuclei: ^^Sr, 
^°Zr and ^^Mo. Calculations within the RQRPA are shown by the dashed curves, and the 
RQTBA - by the solid curves. Experimental data are taken from the EXFOR database 621]. 
These three figures clearly demonstrate how the two-quasiparticle states, which are responsi- 
ble for the spectrum of the RQRPA excitations, are fragmented through the coupling to the 
collective vibrational states. The effect of the particle-vibration coupling on the low-lying 
dipole strength below and around the neutron threshold within the presented approach is 
shown in the left panels of the Figs. [T]121 Our calculations for the tin chain give us an 
example how the low-lying strength develops with the increase of the neutron excess. In the 
doubly-magic ^°°Sn two first relatively weak RRPA peaks appear between 9 and 10 MeV. 
Quasiparticle-phonon coupling redistributes these structures and shifts them about one MeV 
lower. In the ^°^Sn due to the pairing correlations in the neutron system the whole RQRPA 
picture is shifted towards higher energies, and there is practically no strength below 10 MeV. 
In the corresponding figure we find only the strength caused by the fragmentation of the 
higher-lying RQRPA peaks above 11 MeV. In the ^^'^Sn the neutron excess becomes enough 
to form the pronounced pygmy mode situated in the RQRPA at about 9.2 MeV and spread 
over many states of the 2q phonon nature beginning from 5 MeV. Fig. [2] shows how this 
tendency develops in the more neutron-rich nuclei: more strength is split to this region and 
this strength goes to lower energies. 

The Lorentz fit parameters for the calculated GDR in the energy intervals: (10-22.5) MeV 
for the tin chain and (10-25) MeV for the = 50 chain are displayed in Table [T] and they 
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are compared with the corresponding data of Refs. [63|, |6J] . In our work the Lorentz fit is 
performed in such a way that the obtained Lorentzian has the same momenta of -2,-1 and 
zero orders as our microscopical strength function. This method works well if the model 
strength function is rather close to the Lorentz shape. From the Tabled we notice that the 
inclusion of the particle-phonon coupling in the RQTBA calculation induces a pronounced 
fragmentation of the photo absorption cross sections, and brings the mean energies and 
widths of the GDR in much better agreement with the data, for all the examined nuclei. 
The contribution of the low-lying strength below 10 MeV to the dipole spectrum is quantified 
in Table HTl For the each nucleus, we have calculated the following quantities: the non-energy 
weighted sum ^5(£'l) f, which is obtained by direct integration of the strength, and the 
energy-weighted quantity 'Y^ EB{E1) f, which is an integral of the cross section expressed 
in the percentage of the classical Thomas-Reiche-Kuhn sum rule. Both quantities have been 
calculated with RQRPA and RQTBA to emphasize the effect of the quasiparticle-phonon 
coupling for the two energy intervals: (0-10) MeV and (0-8) MeV. The choice of the intervals 
is determined by the fact that, from one hand, in our approach we associate the pygmy modes 
with the dipole strength which originates from the first pronounced RQRPA peaks of the 
isoscalar nature, and these peaks are situated in the tin isotopes just below 10 MeV. From 



the other hand, the measurements o 



photon scattering experiments 65 



" the low-lying strength excited in these nuclei in the real- 



66l | are restricted by the energy around 8 MeV because 



at higher energies the sensitivity of these experiments decreases considerably. Therefore, we 
have included the strength, calculated in the both energy intervals, into the Table HTl 
The integral contribution of the low-energy portions calculated within the RQTBA agrees 
very well with the available data that can be seen from Table HTl the inclusion of the coupling 
to phonons noticeably improves the description. Moreover, below 8 MeV in the most of the 
examined nuclei we observe that the quasiparticle-phonon coupling is the only mechanism 
which brings the strength to this region where the pure RQRPA has no solutions at all. We 
have found also general agreement of our results for ii6,i30gj^ ^^^^ ssgj, relatively 



recent studies of the low-lying dipole strength in Refs. |67|, 168| in the Quasiparticle Phonon 
Model (QPM) {29], although our analysis of transition densities leads to somewhat differ- 
ent conclusions (this analysis will be considered in a special publication). In the QPM the 
model space included up to three-phonon configurations built from a basis of QRPA states, 
calculated with the separable multipole-multipole residual interactions with adjustable pa- 
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rameters, that could be a possible source of the above mentioned differences. In these works 
a clear dependence of the PDR strength and centroid energies on the neutron-skin thickness 
was demonstrated. 

Although the integral strength is described rather good within our approach, the level den- 
sities of the obtained low-lying spectra seem to be underestimated as compared to the ex- 
perimental works. In the other words, in our approach the effect of fragmentation of the 
RQRPA excitations due to coupling to phonons is not enough. To obtain a more realistic 
effect, we could, obviously, use the experience of the previous calculations within the non- 



relativistic models. In the Refs. 



65, 



661], the calculations within the QPM model 29|] with 



taking into account one-, two- and, in the Ref. 65||, also three-phonon configurations lead to 
the better description of the PDR fine structure, although some parameters of the residual 
interaction were fitted. We could include at least more vibrational modes into our phonon 
subspace, that will not even require any modification of the model. Another way to enrich 
the spectrum is to take into account ground state correlations of the singular type, according 
to Ref. 



The fragmentation of the resonances, induced by the quasipartic 
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e-phonon coupling, is a 



49| (see also relatively 



very well known result which has been obtained long ago [47 , 
recent calculations of electric dipole excitations in open-shell nuclei including QPC within the 
framework of non-relativistic approaches based on the Skyrme energy functional 
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5l| as 



we 



1 as on the simple semi-phenomenological scheme including the single-particle continuum 
281]). Actually, one finds more or less a similar level of agreement between the available 
experimental data and the theoretical predictions of these approaches. We notice, however, 
that, in general, our self-consistent relativistic approach reproduces the shapes and often the 
mean energies of giant dipole resonances better than the other above mentioned approaches, 
that could be attributed to the more realistic form of the meson-exchange force and to the 
fully consistent calculation scheme. 

From Figs. [2l[3]one can see that the envelopes of the calculated GDR in ii6,i20gj-^ between 10 
and 22 MeV and in ^°Zr, ^^Si between 12 and 25 MeV are rather close to the experimental 
cross sections. The deviations from the smooth Lorentz shape observed in experiments could 
be attributed to some minor drawbacks of our approach and calculation scheme: neglecting of 
the more complicated, than the 2q^ phonon, configurations by the time blocking, discretized 
continuum, restriction of the phonon subspace by the only low-lying modes, and, at last, too 
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FIG. 1: The calculated dipole spectra for the light tin isotopes. Right panels: photo absorption 
cross sections computed with the artificial width 200 keV. Left panels: the low-lying portions of 
the corresponding spectra in terms of the strength function, calculated with 20 keV smearing. 
Calculations within the RQRPA are shown by the dashed curves, and the RQTBA - by the solid 
curves. [62 1. 
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FIG. 2: The same as in Fig. [H but for heavier tin isotopes, compared to data of Ref. [62] for 
ii6,i20sn. 

simple model for the pairing force. Nevertheless, we find that the agreement with the data 
for the GDR cross sections in these nuclei is very good, especially taking into account the 
fact, that our approach is fully consistent and contains no any fit additionally to the fit of 
the RMF energy functional parameters NL3 which are fixed in the very beginning and used 



36 



RQRPA 
RQTBA 



> 

0) 



1000- 



800 



> 
0) 



RQRPA 
RQTBA 



> 



0) 



PDR ''Mo 



5 6 7 



4 f-i^ 



i 



800- 



600- 



400- 



200- 



RQRPA 
RQTBA 




g 

"-4— » 

o 

w 

w 
w 
o 

o 



c 
o 

"-4—' 

O 


w 
w 

CO 

o 
o 



10 15 20 25 30 





RQRPA 






RQTBA 






GDR ^'Mo 











c 
o 

o 

w 

W 
O 



9 10 11 12 

E [MeV] 



10 



15 20 25 

E [MeV] 



30 



FIG. 3: The same as in Fig. [H but for = 50 isotopes, compared to data from Ref. 64] 



for the entire nuclear chart. Therefore, we conclude that the main mechanisms which are 
responsible for the damping of the GDR are taken into account correctly and consistently. 
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TABLE I: Characteristics of the isovector dipole spectrum for the examined A'^ = 50 and Z = 50 
nuclei: mean energies {E), widths T and EWSR values calculated with the RQRPA and with the 
RQRPA extended by the particle-phonon coupling (RQTBA), compared to data. The values of (E) 
and r have been obtained by a Lorentz fit of the computed strength functions within the interval 
[Sn,<iSn] where 5„ is the neutron separation energy. 
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V. OUTLOOK AND CONCLUSIONS 

The Relativistic Quasiparticle Time Blocking Approximation (RQTBA) has been developed 
and applied for nuclear structure calculations. The physical content of this approach is the 
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TABLE II: Integral characteristics of the isovector dipole spectrum for the examined = 50 and 
Z = 50 nuclei: the integrated strength of the low-lying part below 10 MeV, calculated with the 
RQRPA and with the RQRPA extended by the particle-phonon coupling (RQTBA), compared to 
the available data. 
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quasiparticlc-vibration coupling model based on the relativistic energy density functional 
and the relativistic QRPA. The approach is formulated for a system with an even particle 
number in terms of the Bethe-Salpeter equation for the p/i-channel in the doubled space to 
describe a response of the system in an external field and its spectral characteristics. 
The static part of the single-quasiparticle self-energy is determined by the relativistic energy 
functional with the parameter set NL3 based on a one- meson exchange interaction with a non- 
linear self-coupling between the mesons. An independent phenomenologically parameterized 
term is introduced into the relativistic energy functional to describe pairing correlations 
which are considered to be a non-relativistic effect and treated in terms of Bogoliubov's 
quasiparticles and, in the application, within the BCS approximation. In order to take the 
QPC into account in a consistent way, we have first calculated the amplitudes of this coupling 
within the self-consistent RQRPA with the static interaction. Then, the calculated QPC 
energy- dependent self-energy was introduced into the Dirac-Hartree-Bogoliubov equation for 
the single-quasiparticle wave function and into the equivalent Dyson equation for a single- 
quasiparticle Green's function. The Bethe-Salpeter equation for the response function in 
the doubled space contains the energy-dependent induced interaction connected with the 
energy-dependent self-energy by the consistency condition. The BSE has been formulated 
and solved in both Dirac-Hartree-BCS and momentum-channel representations. 
In order to solve the BSE in the quasiparticle time blocking approximation, the time- 
projection technique is used to block the two-quasiparticle propagation through the states 
which have more complicated structure than 2q phonon. The nuclear response is then 
explicitly calculated on the 2q (8) phonon level by summation of infinite series of Feynman's 
diagrams. In order to avoid double counting of the QPC effects a proper subtraction of the 
static QPC contribution has been performed. Since the parameters of density functional for 
the static RHB description have been adjusted to experiment they include already essential 
ground state correlations. 

The RQTBA introduced in the first sections of this work is applied for the calculation of 
spectroscopic characteristics of the isovector dipole excitations in the wide energy range up 
to 30 MeV for spherical open-shell nuclei, in particular for the isotopes I00'i06,ii4,ii6,i20,i30gj^ 
(Z=50) and ^«Sr, 9°Zr, ^^Mo (N=50). The QPC leads to a significant spreading width of the 
GDR as compared to RQRPA calculations and causes the strong fragmentation of the pygmy 
dipole mode and its spreading to lower energies. This is in an agreement with experimental 
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data as well as with the results obtained within the non-relativistic approaches. 
The good agreement of our results with the experimental data obtained without any addi- 
tional adjustable parameters for a large number of semi-magic nuclei, confirms the univer- 
sality of the RMF energy functional and the predictive power of our approach. We hope, 
that some of the minor drawbacks in these calculations can be overcome by using in future 
an improved version of the density functional both in the ph- and in the pp-channel. 
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